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We address the question of weak versus strong universality scenarios for the random-bond Ising 
model in two dimensions. A finite-size scaling theory is proposed, which explicitly incorporates 
InL corrections (L is the linear finite size of the system) to the temperature derivative of the cor- 
relation length. The predictions are tested by considering long, finite-width strips of Ising spins 
with randomly distributed ferromagnetic couplings, along which free energy, spin-spin correlation 
functions and specific heats are calculated by transfer- matrix methods. The ratio y/u is calculated 
and has the same value as in the pure case; consequently conformal invariance predictions remain 
valid for this type of disorder. Semilogarithmic plots of correlation functions against distance yield 
average correlation lengths £ av , whose size dependence agrees very well with the proposed theory. 
We also examine the size dependence of the specific heat, which clearly suggests a divergency in the 
thermodynamic limit. Thus our data consistently favour the Dotsenko-Shalaev picture of logarith- 
mic corrections (enhancements) to pure system singularities, as opposed to the weak universality 
scenario. 
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I. INTRODUCTION 

In the study of random magnetic systems, a 
frequently-asked question is whether or not quenched 
disorder destroys a sharp phase transition and, in the 
latter case, whether critical exponents are the same as 
for the .corresponding pure magnetsaci The Harris 
criterioncl provides useful guidance in a number of cases: 
if the exponent a, characterising the divergence of the 
specific heat at the critical point of the pure system, is 
positive then randomness induces crossover to a different 
universality class; for negative a the disordered system 
is expected to exhibit the same critical behaviour as the 
unperturbed one. However, such a rule is inconclusive 
for the subject of this work, the two-dimensional Ising 
model, where the specific heat of the pure system di- 
verges logarithmically (that is, with a = ) at the criti- 
cal point. Further, the Harris approach is perturbative in 
the sense that only weak randomness is considered. Non- 
perturbative methods are thus required, especially when 
one wishes to investigate strongly disordered systems. A 
suitable way to deal with this sort of problems is through 
numerical calculations on finite systems. One then has 
to account for finite-size effects before extrapolating to 
the thermodynamic limit. This is done by testing spe- 
cific hypotheses bearing upon the nature of asymptotic 
behaviour. 

In the present paper we investigate the theoretical pre- 
diction (see Refs. ||,||and references therein) that disorder 
affects the phase transition of the two-dimensional Ising 
model only via a specific, well-defined set of logarithmic 
corrections to pure-system critical behaviour; here we ex- 
tend, and give further details of, the results preliminarily 
reported in Ref. [5]. Such a prediction is in contrast to re- 



cent workO~13, according to which critical quantities such 
as the zero-field susceptibility and correlation length dis- 
play power-law singularities, with the corresponding ex- 
ponents 7 and v changing continuously with disorder so 
that the ratio 7/1/ is kept constant at the pure system's 
value (the so-called weak universality scenarioo) . 

Here we calculate free energies and spin-spin cor- 
relation functions on long, finite-width strips of two- 
dimensional disordered Ising systems. The main motiva- 
tion for the use of this geometry is that strip calculations, 
together with finite-size scaling (FSS) conceptsEM are 
among the most accurate techniques to extract critical 
points and .exponents for non-random low-dimensional 
systemsHLj. The rate of decay of correlation functions 
determines correlation lengths along the strip. These 
latter are, in turn, an essential piece of Nightingale's 
phenomenological renormalisation schemeE3E3, and have 
been given further relevance via the connection with crit- 
ical exponents provided by conformal invariance con- 
cepts!-! Early extensions of strip scaling lo r.a jidom sys- 
temsEj have since been pursued furthciij'tjO and put 
into a broader perspective. In particular, it has been 
shown that although in-sample fluctuations of correla- 
tion functions do not die out as strip length is increased, 
averaged values converge satisfactorilyta; throughout the 
present paper we shall make use of this fact to calculate 
error bars of related quantities. 

We consider the two-dimensional Ising model on a 
square lattice with bond randomness. The particular 
version of disorder studied in this work is a binary dis- 
tribution of ferromagnetic interaction strengths for both 
vertical and horizontal bonds, 
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P{Jij) = ^Wa - Jo) + S(Jij - rJ )) , < r < 1 , 

(1) 

which is the prototypical random-bond Ising system, and 
exhibits the unique advantage that its .critical tempera- 
ture (3 C = 1/ksTc is exactly knownllSEl as a function of 
r through: 



sinh(2/3 c J )sinh(2/3 c rJ ) = 1 



(2) 



For given r one can then sit at T — T c {r) and be sure 
that numerical errors due to imprecise knowledge of the 
critical point are absent. Also, a vast amount of simula- 
tional work has been done on this same model (see Ref. 
|), thus comparison is made easier when appropriate. 

The layout of the paper is a*, follows. We first recall 
m Sec. p the main predictionsau concerning logarithmic 
corrections to the singular behaviour of bulk quantities 
for disordered two-dimensional Ising model, and discuss 
how such corrections should show up in the correspond- 
ing finite-size quantities. In particular, we show that a 
logarithmic term is expected to be the leading correction 
to the finite-size be hav iour of the temperature deriva- 
tive of In Sec. Ill we outline numerical aspects of 
our calculational approach for the magnetic susceptibil- 
ity, the correlation length, and the specific heat; also, 
the numerical results are presented and discussed. Sec. 
[V summarizes our findings. 



II. LOGARITHMIC CORRECTIONS AND 
FINITE-SIZE SCALING 

For infinite-system quantities close to the critical point, 
with t = (T — T c )/T c , the following forms have been 
proposed (see Refs. 0,0 and references therein) for the 
correlation length ^ (t) and initial susceptibility Xoo (t) ■ 



r v [i + cin(iA)] 5 , 



Xoo~^ 7 [l + Cln(l/t)]' 



(3) 



(4) 



where v = 1, P = 1/2, 7 = 7/4, 7 = 7/8 and C is 
a disorder-dependent positive constant; for C — one 
recovers pure-system behaviour. Corresponding expres- 
sions have been derived for magnetisation and specific 
heat, which will not concern us for now. Thus theory 
predicts that the dominant power-law singularities (with 
the same indices as for the uniform system) will actually 
be enhanced by logarithmic, divergences. We shall keep 
to current use in the fieldBti and refer to these latter as 
corrections, though strictly speaking the term is inappro- 
priate. 

In searching for signatures of such diverging logarith- 
mic corrections in systems of finite size, one must be care- 
ful about applying recipes used when the bulk singularity 



is purely of a power-law nature. For instance, a naive ap- 
plication, to Eq. (H), of the usual shortcutO t — > to 
extract the size dependence at criticality, would yield a 
correlation length growing faster than L, which clearly 
cannot be true. Instead, one must consider the relation- 
ship between bulk quantities predicted by theory and ex- 
emplified by Eqs. (|^) and (|J), namely 

Xoo ~ (Hoo) 1 ^ ■ (5) 

To see what this implies, recall the FSS hyphotesis00 
for a generic quantity Qz,(f) : 



Q L (t) = f(L) g(z) 



„_ €°o(*) 



(0) 



where L is the linear lattice size and one assumes small 
t, large L. As is well known, the L-dependence must be 
removed as z — > 0. It is immediate that, whenever the 
relationship between Qtx>(i) and £00 (i) is a power law 
such as in Eq. (^|) above, /(£) will be a power law as 
well. This, together with the complementary condition 
that only the L-dependence must remain asymptotically 
for z ^> 1, ensures in the case that the finite-size suscep- 
tibility at the critical point must be 



Xl(0)~L^ = j L 7 / 4 



(7) 



In other words, FSS implies that logarithmic corrections 
must not show up, and the finite-size susceptibility at 
T c will exhibit pure power-law behaviour against L, with 
the same power as in the homogeneous case. The same 
argument is, of course, valid for £l(£) which must then 
scale linearly with L at T c . As shown below, numerical 
data bear out such predictions, for both xl(0) and £l(0). 

This raises the question of how, on finite systems, to 
seek evidence for effects of the bulk corrections predicted 
in Eqs. (||)-(Q). In the following we show that the proper 
quantities to consider are temperature derivatives of e.g. 
£l(£). We apply standard FSS concepts to show that, 
although the dominant behaviour of such quantities is in 
powers of L, the leading corrections to FSS must depend 
on InL. This is in contrast with the corresponding (non- 
diverging) corrections to FSS for, say £ which can be 
fitted by inverse power laws (see Ref. [lj] and below). 

First we recall that the FSS form for £ is, from Eq. ra): 



z< 1 



const., z 1 
The temperature derivative of £z, is then 

WW = — — = (z) , ( z>1 



(8) 



dt 



(9) 



where (ioo(t) = d^oo/dt and the prime denotes a deriva- 
tive with respect to z. While the bulk limit z -C 1 of Eq. 
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(0) is a straightforward identity, the vanishing of 4>'(z) 
for 2 > 1 [as implied by Eq. (||)] must be qualified. In- 
deed, HL(t) does not diverge in the latter regime, while 
Moo(i) does when t — > 0. Thus </>'(.z) ~ (^oo(U)~ , in the 
sense that the dependence of <f>' on t through ^ must 
be such as to cancel the diverging i-dependence of n^. 
Since the FSS ansatz predicts that t only arises through 
the ratio £oo(t)/L, one can deduce the L-dependence of 
Hhit) for z 1. Up to now, the argument is entirely 
general and variations of it have been commonly used in 
the FSS literature. 

Turning to the two-dimensional random-bond Ising 
model, where the bulk quantities are expected to be- 
have as in Eqs. (j^-Q), one has for t <C 1 (consistent with 
our goal of deriving expressions suitable for the z 3> 1 
regime) : 



&o~t-" (ml/if 



(10) 



which can be iteratively inverted to give t as a function 
of Coo: 



The expression for /ioo(i) is 

M*) ~ t~ il+u) [l + Chi(l/t)f, 
plus less- divergent terms, which for f«l can be put as 

/ioo(i) ~ (InlAf 



(11) 



(12) 



W*^ 1 '" (In foe) 



— I// f 



(13) 



where Eq. (11) was used in the last step. It follows im- 
mediately that 



which, when plugged back into Eq. 
©, gives: 



(14) 

together with Eq. 



z > 1, 



(15) 



so all diverging factors related to £oo are removed, but 
a non-diverging ^-dependent term remains which even- 
tually vanishes. Strictly speaking, Eq. (|l5|) means that 
for both t <C 1 and £ » 1, but such that z >■ 1, 
one must observe essentially the leading power-law form 
/il ~ L 1+1 / v . However, even though Eq. (|J) enables 
one to sit exactly at t = 0, Eq. ( |l5| ) suggests the ex- 
istence of a regime in which the leading correction to 
power-law behaviour is ~ (1 — A In L) v l v for finite and not 
very large strip widths L, thus defining an effective (non- 
diverging) screening length £ s = e x l A . This heuristic 
procedure draws on ideas used to interpret experimental 
data for systems where a full divergence of the correLa=j 
tion length isJaindered by percolatioro, random fieldt£l 
or frustrationEj effects. Defining the inverse correlation 
lengths k (actually observed), n a ~ t" and k s = (£ s ) _1 



(representing the physical factor that smears the diver- 
gence, e.g. domain size), one writes 



k = Kq + K s 



(16) 



with good resultsE3t3. Here, n s does not originate from a 
physical feature of the infinite system; instead, it reflects 
the overall effect of higher-order corrections in such pre- 
asymptotic region (strip widths L < 15). While £ s is 
of a different nature to the crossover length _Lp setting 
the scale above which disorder effects are felttm, the two 
lengths vary similarly with disorder, as explained below. 

We now describe the numerical procedures used to test 
the predictions given by Eqs. (0) and ([l5]), and the re- 
spective results. 



III. CALCULATIONAL METHOD AND RESULTS 

We have used long strips of a square lattice, of width 
4 < L < 14 sites with periodic boundary conditions. 
In order to provide samples that are sufficiently repce-. 
sentative of disorder, we iterated the transfer matrixt2l 
typically along 10 7 lattice spacings. 

At each step, the respective vertical and horizontal 
bonds between first-neighbour spins i and j were drawn 
from the probability distribution Eq. ([!]) above. We have 
mainly used three values of r in calculations: r = 0.5, 
0.25 and 0.1; the two smallest values have been chosen 
for the purpose of comparing with recept Monte-Carlo 
simulations where v and jy are evaluatedcil. The critical 
temperatures, from Eq. (g), are: T c (0.5)/ Jo = 1-641 . . . ; 
T c (0.25)/ J = 1.239 . . . ; T c (0.1)/ J = 0.9059 ... (to be 
compared with T c (1)/Jq = 2.269... ). We also evalu- 
ated critical correlation lengths and their derivatives for 
r = 0.01 and 0.001, with respective critical temperatures 
T c j J = 0.5089 ... and 0.3426 .... 



A. Susceptibility 

The calculation of finite-size susceptibility data and 
their extrapolation goes as follows. First, we include 
a uniform longitudinal field h in the Hamiltonian, and 
obtain the largest Lyapunov exponent A° for..ajStrip of 
width L and length N ^> 1 in the usual wayESEZl. Start- 
ing from an arbitrary initial vector Vo, one generates the 
transfer matrices T t that connect columns i and i + 1, 
drawing bonds from the distribution Eq. (Q), and applies 
them successively, to obtain: 



A 



N 



In 



v 



(17) 



The average free energy per site is then fj^ ve (T, h) = 
— -^A^, in units of fc^T. The initial susceptibility of a 
strip, xl(T c ), is given by 
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Xl(T c ) 



d 2 fr e (T,h) 



dh 2 



= IP'" Q(0) , (18) 



T=T c ,ft=0 



where, according to the discussion in the preceding Sec- 
tion, we assume a pure power-law dependence on L at 
T = T C . 

Ati.™ e (T, /i) is expected to have a normal distribu- 
tionE3l£3, so will \l- Thus the fluctuations are Gaus- 
sian, and relative errors must die down with sample size 
(strip length) iV as 1/yN. Typical strip lengths varied 
from N = 2 x 10 6 (for r = 0.5) to N = 2 x 10 7 (for 
r = 0.1), which are much longer than those used in Ref. 
|l7j ; they provide estimates for the free energy with an 
accuracy of 0.01%, which is crucial to compute reliable 
numerical derivatives. In order to get rid of start-up ef- 
fects, the first No = 10 5 iterations were discarded. The 
intervals (of external field values, in this case) used in 
obtaining finite differences for the calculation of numer- 
ical derivatives must be strictly controlled, so as not to 
be an important additional source of errors. We have 
managed to minimise these latter effects by using Sh 
typically of order 10 -4 in units of J when calculating 
f L ave (T c ;h = 0, ±5h) for the derivative in Eq. ©. We 
estimated the first Lyapunov exponent at (T = T c , h = 0) 
and (T — T c , h — ±5h) with four different realizations of 
the impurity distribution, each one giving a separate esti- 
mate of the initial susceptibility. From them the average 
Xl(T c ) and the error bars (twice the standard deviation 
among the four overall averages) are taken. 

A succession of estimates, ("y/isjr , for the ratio 7/1/ is 
then obtained from Eq. rtlS) as follows: 



In [xl{T c )/xl-i{T c )\ 
In[L/(i-l)] 



(19) 



The respective error bars follow from those of the corre- 
sponding finite-size susceptibilities. In order to extrapo- 
late this sequence, we refer to early work on the eigen- 
value spectrum of the transfer matrix for pupa systems 
with a marginal operator in the HamiltoniancJ. There, 
it is shown that the critical free energy per site is af- 
fected only by an additive logarithmic term in the co- 
efficient of the leading, L~ 2 -dependent, finite-size, cor- 
rection (proportional to the conformal anomalycJ, c): 
f(L)-f (00) = - (tt / 6L 2 )[c+B (In L)- 3 +...}. Since disor- 
der is expected to be marginally irrelevant in the present 
case, and assuming that the field derivatives commute 
with the L-dependence (at least as dominant terms are 
concerned), we expect a similar picture to hold here. Of 
course, with the imprecisions introduced by randomness 
one can only expect to see the leading power-law depen- 
dence (see, e.g., Ref. [I?] for further illustrations of this 
point). 

Least-squares fits for plots of (jjv) L against l/L 2 
provide the following extrapolations: 7/1/ = 1.748 ± 
0.012, 1.749 ± 0.008, and 1.746 ± 0.013, respectively for 
r = 0.50, 0.25, and 0.10; the latter two estimates agree 
with 1.74 ± 0.03, 1.73 ± 0.05, obtained in Ref. M. 



The overall picture is thus consistent with the predic- 
tion of Eq. (0), that is j/u — 7/4, same as for the pure 
system, for all degrees of disorder. Recalling the Intro- 
duction, this still is not enough to distinguish between 
weak- and strong-universality scenarios, as both coincide 
in their predictions for the ratio of exponents. One needs 
to try and isolate a single exponent, which will be done in 
the next subsection through investigation of correlation 
lengths. 

Taken together with the results of Ref. [l?] where n was 
found to be 1 /4 through an analysis of averaged correla- 
tion lengths, and using the scaling relation 7/1/ = 2 — 77, 
the present analysis of finite-size susceptibilities gives in- 
dependent suppor|_to the view that: (1) the conformal 
invariance relatiorO rj = L/tt^l(T c ) still holds for disor- 
dered systems, provided that an averaged - as opposed to 
typical, see next subsection - correlation length is used; 
and that (2) the appropriate correlation length to be used 
is that coming from the slope of semi-log plots of cor- 
relation functions against distancdUJ. Interestingly, the 
connection with the conformal invariance prediction also 
rules out any explicit diverging logarithmic L-dependence 
on £ L . 



B. Correlation lengths 

The aim of this subsection is to check on the validity 
of Eq. (|l5|), or rather, its predicted consequences in the 
pre-asymptotic region within our reach, t <C 1, L < 15. 

The first difference to the free energy calculation de- 
scribed above is that the correlation functions are_£x=. 
pected to have a distribution close to /o^normacj'Ej 
rather than a nptmal one. This has been thoroughly 
checked recentlyca. Thus self-averaging is not present, 
and fluctuations for a given sample do not die down with 
increasing sample size. However, it has been numerically 
verified that the spread among overall averages (i.e. cen- 
tral estimates) from different samples does shrink (ap- 
proximately as N^ 1 / 2 ) as the samples' size (TV) increases 
(see Fig. 2 of Ref. flq). Accordingly, in what follows 
the error bars quoted arise from fluctuations among four 
central estimates, each obtained from a different impu- 
rity distribution. Similar procedures seem to have been 
followed in Monte-Carlo calculations of correlation func- 
tions in finite (L x L) systemsEil. 

The direct calculation of correlation functions, (o~oo~r) , 
follows the lines of Section 1.4 of Ref. |l3|, with standard 
adaptations for an inhomogeneous systemEU. For fixed 
distances up to R = 100, and for strips with the same 
length as those used for averaging the free energy, the 
correlation functions are averaged over an ensemble of 
10 4 -10 5 different estimates to yield (croffl)- 

The average correlation length, £ av , is defined by 



(a a B ) ~ exp [-R/D , 



(20) 



and is calculated from least-squares fits of straight lines 
to semi-log plots of the average correlation function as a 
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function of distance, in the range 10 < R < 100. And, 
finally, £ a " is in turn averaged over the different realiza- 
tions of impurity distributions. 

Recall that, as explained in Ref. [Tt], the inverse of £ a " 
is not the same as the difference between the two leading 
Lyapunov exponents, which gives the decay of the most 
probable. Of typical (as opposed to averaged) correlation 
functiontlftZIH. It has been predicted^ that typical cor- 
relations in bulk two-dimensional random Ising magnets 
decay as {(JqUr)^ i?~ 1/,4 (lni?) -1 / 8 , while for averaged 



form inspired by Eq. ( p5| ) with v 
value) and v = 1/2, namely 



ones as in Eq. (20) logarithmic corrections are washed 
away, resulting in a simple power-law dependence. For 
strips one could expect, in analogy, with the case of pure 
systems with marginal operatorsca, additive logarithmic 
corrections to the leading L^ 1 behaviour of typical cor- 
relations: Al - A° L = {w/L)[r) + D^nL)- 1 + . . .} with 
7] =1/4. 

It has been conjectured that the averaged correla- 
tion functions at criticality of the random-.bpnd Ising 
model are identical to those of the pure case ^j-aumeri- 
cally the two quantities are indeed very closeEj'EJ, while 
most-probable and pure-system correlation functions do 
no fit each other so well, thougi-their L-dependence is 
similarly. Given the exact resultH that, for strips of pure 
Ising spins the corrections to the leading L^ 1 behaviour 
of (£ av )~ 1 as given by Eq. (^p| ) depend on L~ 2 , it seems 
reasonable to expect L~ x {i.e. faster than inverse loga- 
rithmic) terms also in the present case. This has been 
shown to work well, with the same x = 2, in Ref. 0. 

We now proceed to testing Eq. (p"5|). We calculate //£ 
at T c [see Eq. (JsJ) ] numerically, from values of evalu- 
ated at T c ± ST, with ST/T C = 10~ 3 . 

Assuming a simple power-law divergence ^ ~ t~ v ' , 
- i.e., ignoring, for the time being, less-divergent terms 
such as logarithmic corrections - we obtain the estimates 
for systems of sizes L and L — 1 : 



1 



In (M-l/Ml-i) 



T=T C 



ln(L/L - 1) 



- 1 



(21) 



This is-slightly different from the usual fixed-point cal- 
culationEj, and is more convenient in the present case 
where the exact critical temperature is known. Our data 
for each pair of (L, L—l) strips have appeared in Ref. [|, 
and we quote here, for completeness, just the extrap- 
olated (against 1/L 2 ) values: v = 1.032 ± 0.031 (for 
r = 0.5; here we have extended the previous calcula- 
tions up to L = 14), v = 1.083 ± 0.014 (r = 0.25), and 
v = 1.14 ± 0.06 (r = 0.10). Taken at face value, these 
data show a systematic trend towards values of v slightly 
larger than the pure-system value of 1, though the vari- 
ation is smaller than that shown in Ref. 21. 



Before accepting this trend as an indication of the 
weak-universality scenario, we must test for corrections 
to pure-system behaviour caused by less-divergent terms, 
as being responsible for the apparent change of v with 
disorder. We then try to check whether our data fit a 



L 2 



(1 - AhxL) 



1/2 



1 (the pure-system 



(22) 



Prior to displaying our results, we recall that the influ 
ence of randomness is expected to show on scales largi 
than a disorder-dependent characteristic length Lca\ 
For L < Lc one should have apparent pure-system be- 
haviour. 

A plot of (/Xi/i 2 ) as a function of lni for differ- 
ent values of r, including r = 1, is shown in Figure ^. 
The pure-system behaviour consists in a monotonic ap- 
proach to a horizontal line, with ever-decreasing slope. 
For r = 0.50 and 0.25 we can see the pure-system trend 
for small L, followed by a clearly marked crossover to- 
wards a form consistent with Eq. (|22j). In each case, log- 
corrected behaviour sets in for suitably large L, exactly 
in the manner predicted by theory: the data stabilize 
onto a straight line with negative slope only for L > Lc, 
which decreases with increasing disorderaH. One may as- 
sume, admittedly with some arbitrariness, Lc for each r 
to be approximately the location of the maximum of the 
respective curve in Figure ^J. This gives Lc — 8, 5, and 
2 respectively for r = 0.50, 0.25 and 0.10 (for r = 0.10 
data for L = 2 and 3, not shown in the Figure, were used 
as well). 

An order-of-magnitude guide to the size of the pre- 
asymptotic region where Eq. (|2|) is expected to hold, 
such that for larger L the pure power-law behaviour pre- 
dicted by Eq. (03) at t = takes over, is the "screening 
length" £ s = e^of Eq. (||). For r = 0.50, 0.25 and 0.10 
one has the approximate values £ s <~ 4 x 10 16 , 7 x 10 4 and 
4 x 10 2 respectively. Though any of these is far beyond 
the largest strip width within reach of calculations, the 
trend against disorder is clearly similar to that of Lc- 

We thus tried stronger disorder (smaller values of r) , in 
order to look for a signature of pure power-law behaviour 
at a feasible L < 15. In Figure || we show (^/L 2 ) as a 
function of InL for r = 0.01 and 0.001. Proximity to the 
percolation threshold is reflected in the large error bars, 
which render our central estimates virtually meaningless 
for L < 5; for larger L, fluctuations are reduced, owing 
to the exponential growth in the number of intra-column 
configurations, so we still can manage reasonable fits in 
that range. Unfortunately, no clear sign can be seen of a 
trend towards a horizontal line. We believe that a con- 
junction of (i) smaller r, (ii) larger L and (Hi) longer 
strip length N would eventually unearth the expected 
stabilization, though we do not feel secure to venture nu- 
merical guesses at this point. 

The above correlation length analysis thus provides us 
with an interpretation of the numerical—data which, it 
should be stressed, is backed by theoryau, without re- 
sorting to disorder-dependent exponents. Nevertheless, 
we have found that the general statistical quality of the 
data does not allow one to distinguish clearly in favour of 
either possibility, in terms, e.g., of least-squares fits. We 
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therefore seek complementary quantitative information 
through the analysis of specific heat data. 



C. Specific heats 



The same theorya'El that gives rise to Eqs. (|3j-|4|) predicts 
that the singular part of the bulk specific heat per particle 
for the disordered Ising model, near the critical point, is 
given by 



C M (t)~(l/Co)ln(l + C ln(l/t)) 



(23) 



where again Co is proportional to the strength of disor- 
der, and the pure-system simple logarithmic divergence 
is recovered as Co — > 0. For Co ^ and f<la double- 
logarithmic singularity arises, whose amplitude Eq. (|23] ) 
predicts to decrease as disorder increases. The bulk spe- 
cific heat cannot then be put as a simple function of 
the correlation length given in Eq. (|3|), and one cannot 
predict pure-system behaviour against L for finite sys- 
tems, as was the case for the susceptibUity and correla- 
tion length above. Instead, theory givesB 



C L (t = 0) ~ Ci + oln(l + bhiL) 



(24) 



where, similarly to Eq. (]23j), b — > for vanishing disorder. 
In this latter limit the product ab must remain finite, but 
it is not a priori obvious from theory whether a oc b~ x 
away from that. In fact, the form Eq. (24) has been veri- 
fied by Monte-Carlo simulations of L x L systemsa with 
the result that the slope of plots of Cl against lnlnL 
decreases for increasing disorder. This shows that in this 
case the simple FSS recipe t — > L _1 seems to work satis- 
factorily. 

An investigation of the specific heat on strips is clearly 
of interest, in order to check the consistency of our own 
correlation-length data, and also to provide a compari- 
son with the trends found for the specific heat in L x L 
systems, both as described above and in recent worktl 
where a non-diverging behaviour is apparently found in 
the thermodynamic limit. 

Our results are displayed in Figure |[ where one can see 
that the fit to a double-logarithmic form is reasonable; 
for small disorder r = 0.5 we get an overall better fit 
to a pure logarithmic divergence, similarly to the result 
for L x L lattices". This is again because, as disorder 
decreases one gets apparent pure-system behaviour for 
relatively large L. 

The slope of the plots turns smaller for higher disor- 
der, again in agreement with the trend found for L x L 
latticesEj; however, no sign of an eventual trend towards 
non-divergenceu can be distinguished. 

The recent claims that for strongly disordered Ising 
systems in two dimensions, the specific heat is finite at 
T c , have been made on the basis of numerical simula- 
tions of site-diluted modelsB. Specific heats were plotted 
against t (t > 0) for system sizes and temperatures such 



that L/£oo(t) > 1 (thus excluding the very close vicin- 
ity of the transition); see Fig. 1 of Ref. ||. While for 
impurity concentration c = 1/9 a divergence was clearly 
seen, data for c = 1/4 and 1/3 were interpreted as sig- 
nalling a finite bulk specific pheat at the transition. Such 
findings have been criticisecH. At this point it is worth 
recalling experimental data. First, in bulk systems the 
specific heat exhibits a broad regular background against 
which the singular part must be singled out. Early ex- 
periments on the two-dimensional site-diluted Ising sys- 
tem Rb2Co x Mgi_2;F4 showed that the amplitude of the 
singular part of the magnetic specific heat decreases as 
dilution 1 — x increasesEH. However, owing to experimen- 
tal difficulties, chief among them the smearing of T c due 
to sample inhomogeneities, clear peaks could be found 
only for 1 — x < 0.11. Later, results from the more accu- 
rate technique of birefringenceEM-and with a presumably 
higher-quality sample confirmedE3 that the specific heat 
diverges for 1 — x = 0.15, apparently with a single loga- 
rithmic dependence identical to that of the pure system; 
this was ascribed to an extreme narrowness of the region 
where disordered (double-logarithmicXpbehaviour would 
show up (in agreeement with theoryoE3). Though, to 
our knowledge, a systematic study of the variation of spe- 
cific heat of two-dimensional Ising systems against dilu- 
tion, by e.g. birefringence techniques, has not been done, 
useful hints may be taken from the corresponding three- 
dimensionaL-ease of Fe^Zni-a^. There, birefringence 
experimental!] show that as dilution increases, the rela- 
tive position of the (narrow) peak at T c against that of 
the (broad) maximum of a short-range order background 
contribution switches from higher to lower temperatures. 
This fact is not directly related to the particular three- 
dimensional features which are used to explain-the dilu- 
tion dependence of the specific heat amplitudecJ. Thus, 
it is not unlikely that for two dimensions too the appar- 
ent non-diverging behaviour seen, for T > T c at c = 1/4 
and 1/3, represents only the background. To see the ac- 
tual (probably small) peak one would have to go closer 
to T c ; imprecisions in the knowledge of T c for site-diluted 
systems (see e.g. Ref. |l6]) may be of capital importance 
then. 

In contrast to this, here and in Ref. |34| one sits right 
at the exactly known T c . Further, according to the dis- 
cussion of finite-size specific heats above, the amplitude 
of the peak at the bulk transition translates directly into 
the slope of the plot of Cl against lnlnL, so the regular 
background is easily dealt with. 

In short, the evidence presented here clearly indicates 
that the specific heat diverges at the transition, with 
a double-logarithmic behaviour. Thus the critical ex- 
ponent a is non-negative. Through hyperscaling argu- 
ments, this ties in with our findings for the correlation 
length, as shown in the following. For weak enough dis- 
order, there should be no question about the dimension- 
ality of the system, as opposed to near the percolation 
threshold in the diluted case [corresponding to r = for 
the bond distribution Eq. (|])], where one might argue in 
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favour of substituting the fractal dimension for the actual 
lattice dimensionality. Therefore, hyperscaling should be 
fully applicable with d = 2, which yields 

2 = 1 - ( 25 ) 

Since our specific heat data implies a > (most likely 
a = 0), one must have f < 1, thus excluding the disorder- 
varyins gKBonents given in Sec. Illb and in previous 
works .ETOEd 



IV. CONCLUSIONS 

We have addressed the question of strong versus weak 
universality in the two-dimensional random-bond (i.e., 
exchange couplings being either J or rJ with equal prob- 
ability; < r < 1 measures the degree of disorder) Ising 
model, through extensive transfer matrix calculations. A 
key ingredient in the analysis of our data has been the 
consideration of subtle finite-size scaling (FSS) effects; 
these come about as a result of constraints imposed by 
the Dotsenko-Shalaev theorya for logarithmic corrections 
in the thermodynamic limit. We have established that 
while the correlation length (and the susceptibility) itself 
should display no signature of size-dependent logarith- 
mic corrections, its temperature derivative, [ir, = d^/dt, 
shows a In L dependence (L is the strip width) over a 
wide range of system sizes. Actually, at the (exactly 
known) critical temperature for the infinite system and 
for constant disorder (i.e., fixed r) the behaviour with 
linear size is as follows. For L < Lc, with Lc be- 
ing a crossover length, the system behaves as in the 
pure case; Lc decreases monotonically with disorder and 
2 < Lc < 8 for the values of r we considered. Above 
Lei Ml is dominated by a InL enhancement over the 
usual pure system power law; that is, the numerical data 
can be explained through consistent theories, without re- 
sorting to disorder-varying critical exponents. The FSS 
theory developed here also suggests that as L increases, 
beyond a (heuristically introduced) screening length £ s , 
one will eventually reach an asymptotic regime where the 
logarithmic enhancements will vanish, leaving only pure 
power-law (pure-system-like) behaviour; see Sec.[fi} This 
coherence length tracks the crossover length, in the sense 
of decreasing with increasing disorder, but its order of 
magnitude is way beyond the reach of our numerical ca- 
pabilities (£ s > 10 2 ) for us to venture a more refined 
analysis of this issue. Note, however, that when t — ► 
after the thermodynamic limit has been taken (which is 
an entirely different matter) it is expected that ln(l/£) 
corrections, as predicted by the Dotsenko-Shalaev the- 
ory, should manifest themselves. Also, our data inde- 
pendently confirm that the conformal invariance result 
£ av = L I ttt] is still valid for the two-dimensional random- 
bond Ising model, with rj = 1/4 as in the pure case. 

As a further test on the consistency of the proposed 
scenario, we have examined the size-dependence of the 



specific heat for this system. Consistently with the above 
findings, the specific heat was seen to be clearly divergent 
in the thermodynamic limit. Since there are no phys- 
ical grounds to invoke a mechanism leading to changes 
in the hyperscaling relation, the case for weak universal- 
ity cannot be supported by our data. Further, it must 
be noted that-a variety of studiea_a£,Lhis problem, both 
theoreticaOo and experimentaOEja concurs with the 
idea that the leading singularities remain the same as in 
the pure case, though they have not dealt with the de- 
tection of logarithmic corrections. 

As regards works whose,|-C|Onclusion is that weak- 
universality holds insteadQ'LlO, though £l(T) and the 
susceptibility xl{T) were calculated, no attempt seems 
to have been made to fit the corresponding data to a 
form similar to Eq. (j22]). Thus it remains to be checked 
whether they would also be consistent with suitable FSS 
expressions based on strong universality concepts. 
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FIG. 1. Finite-size scaling plots of logarithmic corrections 
[Eq. (^)] . Straight lines are least-squares fits of data respec- 
tively for L = 9 - 14 (r = 0.5); 7- 12 (r = 0.25) and 4-12 
(r = 0.1). The error bars are smaller than the data points. 



FIG. 2. Finite-size scaling plots of logarithmic corrections 
[Eq. (0)] for strong disorder. Straight lines are least-squares 
fits of data respectively for L — 6 — 12 (r = 0.001) and 7 — 12 
(r = 0.01) . 



FIG. 3. Specific heat per site at criticality for L = 4 — 12 
and r = 0.50 (squares), 0.25 (crosses) and 0.1 (triangles), 
against mini [Eq. Q]. 
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